فصلنامه فیزیک زمین و فضا
سال چهل و نهم شماره 2 (تابستان 1402)

Studying the bedrock geometry in mining and oil exploration operations to obtain its 2D pattern requires nonlinear reverse computations. Local optimization methods for solving nonlinear inverse problems are based on linearizing the changes of the model similar to a primary model and finding an objective function of minimum error from the model’s parameters; however, these optimization methods are not able to select a suitable primary function that is close enough to the general optimal value. That is to say, every objective function can have several minimum and maximum solutions. The lowest minimum is called the global minimum while the rest of them are named local minimums. Therefore, in local inverse methods, the objective is to find the minimum of an objective function, and also an objective function might have a few local minimums with different values. In this case, it is not suitable to use gradient-based methods for exploration purposes, unless the primary model is very close to the actual answer; which is outside the control of geological structures or the geometry of the subsurface. Despite the easy execution and high convergence rate of the local methods, there is the possibility of being trapped in local minimums because these methods are dependent on the primary model and also finding more than one optimized point in 2D or 3D simulations; this is why local optimization methods are considered deterministic algorithms. Multi-objective and single-objective metaheuristic optimization algorithms are capable of searching the feasible region and they also provide a solution independent of the primary model. Searching the feasible region means finding all the feasible solutions for a problem and each point in this region is representing a solution that can be ranked based on its value. One of the important differences between local optimization and metaheuristic methods is constraining. Constraining metaheuristic global optimization methods are only used for constraining the feasible region based on previous knowledge or estimation relations; which is very different from constraining local optimization that is used for stabilizing inverse simulation. The algorithms used in the present work included a non-dominated sorting genetic algorithm (NSGA-II) and single-objective genetic algorithm, which were used to estimate the depth. The NSGA-II is commonly used to solve problems with multiple, typically conflicting, objective functions. This algorithm is capable of being developed and also has a high potential for solving unbounded multi-objective problems. In addition, the single-objective genetic algorithm (GA) is capable of modeling and solving complex problems. In the present study, both algorithms were verified and validated using the data produced by an imaginary and complex synthetic model. In order for a more precise examination of the performance of both algorithms, the imaginary synthetic data were used both with no noise and with up to 10% Gaussian white noise (GWN). Accordingly, the modeling results indicated a good consistence between the algorithms and the primary model; so that, the root mean square error parameter for the data obtained from the initial data of the synthetic model ranged from 0.05 to 0.35mGal for the NSGA-II and from 0.07 to 0.52mGal for the GA. Also, this parameter didn't exceed 72.4 in the NSGA-II and didn't exceed 93.8 in the GA. Based on the gravimetric modeling of the Western Anatolia, Turkey, the results obtained from both algorithms under similar conditions in terms of parameter settings and number of algorithm executions indicated good performance of the NSGA-II algorithm compared to the single-objective algorithm.

Magnetic Resonance Sounding (MRS), as a surface geophysical method, provides good information about the hydro-geophysical parameters (such as water content and hydraulic conductivity) of aquifers. The main advantage of the MRS method compared to other geophysical methods is that the surface measurement of the MRS signal responds directly to the presence of water below the surface. Despite the high efficiency of this method, the recorded signal is strongly affected by electromagnetic noises, including spike noises and harmonic noises. The first generations of MRS instruments were single-channel instruments. In single-channel instruments, both magnetic resonance excitation and signal recording are done with a single loop, and it is necessary to use various forms of filtering to eliminate noise, in particular powerline harmonics.Then a new generation of multichannel MRS instruments with multiple is was built, so that the main loop is still used for magnetic resonance excitation and signal recording. In addition, a number of reference loops, physically displaced from the main loop, measure only noise. Parts of the noise recorded by the reference loops can correlate with the noise in the main loop. With proper signal processing, the noise in the reference coils can be filtered to obtain a replica of the noise in the main coil and when this replica is subtracted from the signal recorded in the main loop, hence, the MRS signal remains noiseless. One of the current challenges of MRS signal processing is the existence of harmonic noise with the variable fundamental frequency. If the harmonic signal from a specific source has a fundamental frequency that varies with time, most of the proposed algorithms will not perform well in eliminating harmonic noise. Therefore, a new topic that is followed in this paper is to evaluate the performance of the proposed algorithm in cases where the harmonic signal has a fundamental frequency that varies with time.In this paper, in order to obtain an accurate estimate of the parameters of the magnetic resonance sounding signal, a method for eliminating spike and then harmonic noise in the time domain is presented. Synthetic signals are contaminated with different electromagnetic noises to investigate the effect of different optimal filter parameters for spike and harmonic event elimination methods. Spike noise has a detrimental effect on the performance of the harmonic noise elimination algorithm. Hence, spike signals must be deleted or adjusted before applying the harmonic noise elimination algorithm. First, a statistical processing algorithm based on the signal-dependent rank-order mean (SD-ROM) filter for eliminating spike noise is presented and after deleting spike noise, a method for eliminating harmonic noise is assumed based on the fixed and variable fundamental frequency with time using remote reference loop. Numerical results of applying the proposed processing algorithms in the time domain show that by applying the mentioned methods, a considerable amount of spike and harmonic signals are removed leading to a good estimate of the parameters of the magnetic resonance sounding signal (i.e., initial amplitude, relaxation time, phase and frequency of the signal).

Fault slip rate distribution plays an important role in earthquake studies. Faults are loaded at very slow rates in continental interiors. So, interaction among faults and resulting slip distribution can give rise to earthquakes on other faults after a long period of quiescence and seismicity that can migrate from one fault to the onother one.NW Iran-Eastern Turkey is a region of active deformation as a result of oblique collision of Arabia-Eurasia tectonic plates. In northwestof Iran, deformation between the Central Iranian block and the Caucasus domain is accommodated by a fault system and mainly by right lateral strike-slip on the North Tabriz fault. In the current study, we did slip rate partitioning in the fault system of northwest Iranian plateau using the concepts of dislocation theory. Modelling approach is described by Gomberg and Ellis (1994), Flerit et al., (2003) and Armijo et al., (2004) and it differs from rigid block models (Reilinger et al., 2006; Djamour et al., 2011) in which dislocation conditions at the boundaries of blocks are often incompatible with geological evidences. In the alternative method of Flerit et al., (2003), slip everywhere has a direction of motion consistent with geological constraints. The dislocations do not divide the region into closed rigid blocks and slip can vary along strike as observed geologically. Finally we obtain a tectonic model for NW Iran-Eastern Turkey that is more realistic than rigid block model (Reilinger et al., 2006; Djamour et al., 2011) or models based on seismic or geologic strain rates (Haines, 1982; Haines and Holt, 1993; Jackson et al., 1995; Masson et al., 2005). For this purpose we use a three dimensional boundary elements method. First, we consider an elastic and homogeneous half-space for the study area. Then geometric data of fault system are collected from geological and geophysical sources including fault length, width, dip, and locking depth. For Lame coefficients, we use global average values. Both mentioned geometrical and physical data are kept fixed in the modeling process. Then, strain tensor that best fits the GPS data is estimated for the study area using least squares method. Then, stress rate tensor is estimated using generalized Hook’s law. Geomerical chracteristics of faults, physical characteristics of crust and stress rate tensor act as boundary conditions in the model.Faults are locked in normal direction but they are allowed to slip freely in strike and dip directions under the influence of boundary conditions. Regarding the strike changes of faults, the fault surfaces are divided by different segments in strike direction with constant strikes and dips. Then fault segment surfaces are divided into 1km elements. Finally, we have free slipping elements in strike and dip directions as inputs for modeling.Our model is fitted to the fault traces data set of NW Iran-eastern Turkey. The results indicate the dependency of the partitioned slip rate on the boundary conditions and confirm the existence of interaction among faults. Also, partitioned slip rates show that the Chalderan, Guilato-Siahcheshmeh-Khoy, Nakhchivan, North Tabriz and Pambak-Sevan-Sunik faults are right-lateral strike slip in all cases. Also, the slip rate in these faults is almost symmetric and reaches its maximum value around the center of the faults. We show that the maximum value of slip rate in the fault plane is reduced by partitioning, which it will be definitely closer to reality. According to the gridding for slip rate partitioning in the fault system, the highest value of slip rate is always related to the North Tabriz Fault.Previous studies show that the geological slip rate estimates are lower than the present-day GPS-derived slip-rates along the North Tabriz fault. We show that slip rate partitioning solves this discrepancy by considering the mechanical interaction among faults. Our partitioned slip rates for North Tabriz Fault are lower than geodetic rates and are more consistent with geological rates. Finally, we present a model that fits best with the geological constraints.The proximity of the partitioned slip rate to the paleo-seismic values indicates the closeness of the partitioning results to reality with the Boundary Elements Method, compared to other analytical and numerical methods. This research may open new research direction to handle the differene between geologic and geodetic slip rates values in the Iranian Plateau.The boundary elements method is both faster and more accurate for modeling compared to the finite element method used by Khodaverdian et al. (2015). Considering the effect of topography and sphericity of Earth, using the Galerkin boundary element method developed by Thompson (2019) is proposed to get more realistic results. The coefficients matrix in the of Boundary Elements Method is fully populated and in high dimensions it takes a lot of time to solve the resulting system of equations. Sparsing of the coefficient matrix using wavelet transforms is suggested (Ebrahimnejad et al., 2010) in this study. The use of iterative computational methods along with parallel processing will also reduce the computational time (Thompson and Meade, 2019).

The Mohorovičić discontinuity, often known as the Moho, marks the boundary separating the Earth's crust from the mantle. Techniques such as isostatic-gravity and seismic methods can be used to determine this division. The Moho marks the boundary between the continental and oceanic crust and the upper mantle. Simply put, the Moho acts as a physical/chemical boundary between the mantle and the crust and causes significant changes in geophysical properties such as seismic wave velocity, density, pressure, and temperature (Mooney and Masters, 1998; Martinck, 1994; Bagherbandi, 2011 and Dashtbazi et al., 2023). An accurate and high resolution Moho depth model in fields such as geodesy, geology, geophysics, geodynamic modeling, seismic risk assessment, stress field modeling caused by mantle convection (Li et al., 2018; Behr et al., 2022; Singh and Yadav, 2023; Heilman and Becker, 2022; Hashima et al., 2016; Eshagh et al., 2020; Eshagh, 2015 and Gido et al., 2019), and understanding seismic source mechanisms is important, among other applications (Gido et al., 2019 and Dashtbazi et al., 2023). Furthermore, a reliable Moho model can reveal details of crustal structure that provide valuable insights into the complexities of deeper mantle layers; related to the calculations and detailed examination of gravity, geothermal, geomagnetic models (Stalk et al., 2013). Although there are several Moho models, their accuracy and resolution are insufficient in the complex tectonic geometry of the Makran subduction zone (Brizi et al., 2021 and Heilman and Becker, 2022), because these zones show a complex Moho configuration (Shad Manaman et al., 1390; Taghizadeh Farhamand et al., 2015 and Dashtbazi et al., 2023). As a result, the existing models lack the necessary accuracy for the Makran subduction zone, a region approximately 1000 km long located in southeast Iran and southwest Pakistan (Byrne et al., 1992; Shad Manaman et al., 1390, Penney et al., 2017; Dashtbazi et al., 1398; Dashtbazi. et al., 2023).In geophysical and geodetic studies, hybrid methods are mainly used to determine the Moho depth when seismic data with appropriate distribution and abundance are not available. These techniques include the Parker-Oldenberg method and the Wenning-Mines-Moritz method. In an effort to strengthen the existing Moho depth models in the Makran subduction zone, two distinct models named BC and SC through the integration of gravity (VMM) and seismic (CRUST1.0) data, which are processed through the Butterworth filter, spectral combination approaches and the least squares technique, was developed (Bagherbandi , 2011 and Dashtbazi et al., 2023). The resulting models provide a resolution of 5' x 5' degrees of arc, corresponding to a grid size of 9 x 9 km (Dashtbazi et al., 2023). The accuracy of these models was evaluated against four separate regional and local models. The resulting RMS values were 5.28, 1.55, 4.18, and 1.27 km for the BC model and 5.59, 1.17, 3.74, and 3.04 km for the SC model. Also, the Moho depth model obtained for the west Makran region in Iran significantly improved the accuracy and resolution of the Moho depth models in the studied area. The SC Moho model exhibits improved RMS metrics compared to the combined BC model, so we recommend it as the first priority. While the Moho depth models in our research really bring significant improvements to the existing models of the Makran subduction Moho zone, the integration of more detailed seismic data with SC and BC Moho models can improve the developed model for the Makran subduction zone. In the end, we suggest that a similar approach be adopted for the analysis of the Moho model in the eastern Makran region in Pakistan, which allows a comparative evaluation of the Moho depth and structure between the western and eastern parts in order to obtain a better picture of the Moho depth model of the Makran subduction zone.

Fault plane solution is one of the most important tools to determine the orientation of the stress field. The focal mechanism of earthquakes can be applied to determine the direction of rupture propagation, the structure of the fault and the stress field of the region. Dorouneh fault is the largest fault in Iran after the main Zagros fault with about 800 km length. The purpose of this study is to investigate the seismicity, the stress fields and the focal mechanism of earthquakes that have occurred across the three main segments western, middle and eastern parts of the Dorouneh fault. Therefore, the calculation of the focal mechanism is performed using the P-wave first-motion polarity. Also the stress situation of the events and the recognition of fault planes are presented in this research. Along this fault, the blocks have moved in both left and right directions, but certainly one of its last movements has had a right lateral motion. In this study, the waveforms recorded at the seismic stations of the National Seismological Center of Tehran University (IRSC), the National Center of Broadband Seismic Network of Iran belong to the International Institute of Earthquake Engineering and Seismology (IIEES) and seismic stations of Seismological Research Center at the Ferdowsi University of Mashhad were used. At first, the waveforms of each seismic event were combined with each other, and then the relocation of the events were determined based on new data set of this study (Figure 2). A number of high quality earthquake with a magnitude of more than 3 and an average depth of 14 km have been selected (their list is presented in table1) to calculate the fault plane solution.In order to calculate the mechanism of the recent earthquakes in the Khorasan region, interesting results have been obtained based on several methods (Assar Enayati, 1400). One of the common methods for estimating the mechanism of earthquakes, especially earthquakes with small magnitudes and at close distances from the epicenter, is the polarization of the first arrival of the P wave. Due to the dependence of the amplitude and polarization of the P wave on the focal mechanism, by determining the polarization of the first arrival of the seismic phase, the earthquake focal mechanism can be calculated. The results of focal mechanism solutions for significant events around Dorouneh fault show mostly left lateral strike slip motion which is consistent with the tectonic setting of the region. The difference in the focal mechanism of the events in the eastern and western parts of the fault is justified by the northward movement of the Lut block. The integration of data shows high accuracy in calculating the focal mechanism and more certainty about the results. Therefore, in the direction of the Dorouneh fault, the movements are both left-handed and right-handed from the seismic situation. The change in thefocal mechanism obtained from the seismic results, considering the stress axes changes, can represent the second and third order stress fields that balance the stress in this area today. The second-order stress can be related to continental rifting, isostasy adaptation, topography and deglaciation, and the third-order stress field can be related to the local stress source on a scale smaller than 100 km, which is influenced by the structural geometry, and interference between fault systems, topography and local density difference (Sheikh-ul-Islami et al., 2021).

Due to the lack of sufficient information about the density of subsurface masses, in the modeling of the earth gravity field, usually a constant global average density value is used for the entire studied area. However, any accuracy improvement of the mass density used in the modeling of topographic gravitation will increase the accuracy of gravity field modeling. To approve this quantitatively, a topographic density model with a resolution of 30×30 arc second is prepared from the processing of seismographic maps and satellite data of the layers of the lithosphere and used instead of a constant mass density in four study regions inside Iran with diverse topographic and data distribution. These four regions have the dimension of 2.5×3 degrees. The point intervals in the first and third regions (R1 and R3) are approximately 5 minutes, while in the second and fourth regions (R2 and R4) the intervals are approximately 13 minutes. In areas with the same point distribution, R1 has a relatively smoother topography than R3. The topography in R2 is relatively rougher than R4. In addition to the global average density value, the average value for Iran and the region is also included in modeling the gravity topographic masses. Owe to topography diversity, these areas seem to be suitable for an investigation of the RTM technique performance. In modeling the gravity field, the least squares collocation method and, consequently, the RTM technique are used in modeling the effect of topographic mass gravitation. In order to evaluate the effect of lateral density variations when using the RTM technique for gravity field modeling of the earth, the least squares collocation method is used in this research. The RCR technique is used for gravity field modeling by least squares collocation method and the effect of the global topography is removed. To remove the global gravitational effect from the observations, the EIGEN6C4 model degree and order 360 is used. To remove the effect of topography by RTM method, a digital elevation model with a point density of 1 second arc is used. In addition to density correction, the use of improved covariance algorithm in gravity field modeling is also evaluated in this research. The results show that in the areas with more topography, and hence more density variations, the effect of density modification in removing the effect of topography from the gravity anomalies signal of the region is more significant. Furthermore, comparison to the control points of this study shows that application of density correction and the improved covariance algorithm in areas with rough topography and lack of sufficient gravimetric data and proper distribution, the accuracy of gravity field modeling can be improved by 15.6%. Using the IC approach in four regions leads to an increase in modeling accuracy. Among these regions, the highest increase in accuracy is related to region 2 (with a decrease of 1 mGal of standard deviation) and relatively related to region 3 (with a decrease of 11.5% in standard deviation). These two areas have rougher topography than areas 1 and 4.

Shear wave velocity is a key factor to estimate the elastic and petrophysical parameters of the hydrocarbon reservoir. However, shear wave velocity is rarely logged at wells due to the imposition of high costs. Therefore, it is usually attempted to estimate this parameter by different methods from the available and related data. Describing the elastic parameters of reservoir rock, including shear modulus, bulk modulus and Poisson's ratio, requires the measurement of density and compressional and shear wave velocities of the reservoir formations. Direct measurement of the shear wave velocity is done by drilling cores and DSI (Dipole Shear Sonic imager) tools, which are unfortunately very time-consuming and expensive. In this study, a practical method for estimating shear wave velocity in a sandstone oil reservoir is presented. In the studied reservoir, from seven existing wells, the shear wave velocity has been measured by DSI tools in only one of them (well #7). The shear wave velocity log in the location of the other wells was estimated using a petrophysical equation, defined for the location of well #7. The correlation of other logs (i.e. acoustic, density, porosity, resistivity, gamma ray, dolomite volume, quartz volume, and water saturation logs) with the shear wave velocity was investigated in well #7. We found that the compressional wave velocity, density, porosity, dolomite volume and quartz volume logs were more correlated with the shear wave velocity log in well #7. Thus, these logs were selected as input for estimating shear wave velocity log and the experimental equation using the multivariable linear regression method was calculated. The estimated shear wave velocity log using the obtained relationship has a 90% correlation with the measured shear wave velocity log in well #7. Using this petrophysical relationship, the shear wave velocity were estimated in the other wells (blind wells). The main goal in this study, was to produce the volume of the shear wave velocity information at the sandstone reservoir. To obtain 3D volume of shear wave velocity distribution in the reservoir, the seismic and well data are integrated. To achieve this goal, the model-based seismic inversion technique has been performed to obtain the acoustic impedance volume for the sandstone reservoir. The calculated acoustic impedance volume using model-based algorithm has an average of 99% correlation and 15% error with the real acoustic impedance log. The results of the seismic inversion were fed into the cross validation method to derive the optimal number of seismic attributes relevant to shear wave velocity information. The cross validation method shows that the attributes of the filter 20/25-30/45, the cosine instantaneous phase, the acoustic impedance and the instantaneous frequency have the reasonable correlation with the shear wave velocity information respectively, and are selected as the input attributes for the estimation of shear wave velocity volume in the sandstone reservoir. Our results show a good agreement between the real shear velosity log and the predicted shear velocity from the seismic attributes in the place of well #7. The obtained shear wave velocity volume accompanied by the compressional wave velocity information can be used to infer more robust petrophysical parameters in the reservoir.

The potential field data like Gravity data, Magnetic data, Self-potential data, and other natural source data are accessible but hard to interpret and model. Numerous research presents ideas for modeling the potential data; however, they are mostly based on inverse or forward modeling needing a priori constraints and information. In complicated geology and tectonic setting, we do not have convenient access to a priori information to define the constraints. So, we must develop a pure geophysical interpretation method without geological constraints and minimum complexity. In this research, the Normalized Full Gradient’s ability to find the gravity anomaly model was studied. The Normalized Full Gradient method is an effective method for determining anomalous bodies, such as the distribution of oil and gas fields or structural boundaries. The Normalized Full Gradient method depends on the downward analytical continuation of normalized full gradient values of gravity data. Analytical continuation discriminates certain structural anomalies which cannot be distinguished in the observed gravity field. The Normalized Full Gradient of the gravity anomaly is often used rather than the gravity anomaly itself for detecting underground spaces because it is stable and indicates the locations of source bodies. The weakness of the Normalized Full Gradient is the 3D modeling limitation, as we can only calculate the 2D response in practice. On the other hand, the responses can not describe the negative and positive parts of the anomaly. But a unique advantage of the Normalized Full Gradient is, that it does not need the primary information for gravity data modeling. In this research, we used the Normalized Full Gradient for the large-scale Bouguer gravity anomaly interpretation. Bouguer gravity anomaly wavelengths contain information about density distributions of upper mantle and lithosphere structures. A gravity profile is most often a combination of relatively sharp anomalies that must be of shallow origin and very deep and large anomalies with a regional nature.The study profile has a 400 (km) length from SW to NE of North Western Iran, and Sahand, North Tabriz Fault, and Sabalan are the most important structures in the study area. The Normalized Full Gradient synthetic model data study provides the opportunity for the real data recovered model judgment. So, we first showed the Normalized Full Gradient recovered model of the synthetic data test and then based on the resolution of the Normalized Full Gradient used, it provides the lithospheric density of North Western Iran. The result shows the low-density mantle wedge which is probably is beneath the North Tabriz Fault that is responsible for the formation of distinct lithosphere conditions. So, the wedge can explain the complicated tectonic setting of North Western Iran. The mantle wedge has more than 50 (km) wide and more than 40 (km) depth. seemingly, this mantle wedge directly affects the North Tabriz Fault, Sahand, and Sabalan in the shallower depth. The sharpest effect is for the North Tabriz Fault in the shallower part of the mantle wedge and following the shape of the wedge, we can see a sharper effect on the Sabalan in comparison with Sahand.

The boundary layer is the lowest part of the Earth's atmosphere and is directly affected by the Earth's surface. Humans live within the boundary layer, and the height of this layer affects air quality and their physical and mental health. It also affects commercial activities such as air and sea transport. Therefore, measuring the height of the boundary layer and knowing how it evolves is very important.The height of the atmospheric boundary layer is measured in two general ways: in-situ measurements and remote sensing methods. The main disadvantages of in-situ measurements are their cost and the impossibility of continuous measurement. On the other hand, remote sensing methods use the concept of the interaction of sound or electromagnetic waves with the components of the atmosphere. It uses instruments such as Sodar, Radar, and Lidar to continuously study the boundary layer's height at a lower cost.In this study, the data obtained from the elastic backscatter Lidar measurements have been used to extract the boundary layer's height for the city of Zanjan. The Lidar is located at the Institute for Advanced Studies in Basic Sciences (IASBS), and its 532 nm wavelength channel is used in this study. Also, the method used in this study is the Wavelet Covariance Transform (WCT) method. The attenuated backscatter coefficient or the range-corrected signal is expected to decrease for ground-based Lidars as the height increases. As it exits the boundary layer, the concentration of atmospheric particles decreases sharply, and therefore the backscattered signal is expected to drop significantly. The WCT method uses this sharp reduction of the signal to detect the height of the boundary layer. The wavelet covariance transform function, in simple terms, calculates the difference between the Lidar range-corrected signal for a given height interval and the Lidar range-corrected signal for a higher altitude interval. As a result, the peak height equivalent of this function Corresponds to the height of the boundary layer.In this study, only the period 2011-2012 was studied because of the significant number of measurements. In total, we had 105 days of data during the mentioned period, the first of which is related to April 13, 2011, and the last of which is related to October 15, 2012. All measurements are taken for days with no active synoptic conditions. In winter, we have the lowest daily average value of the boundary layer height (0.975±0.556 km), and in summer, we have its maximum value (2.597±0.714 km). In the two seasons of spring and autumn, the extracted values of the height of the boundary layer are very close to each other and are about 1.9 km above ground level (AGL). Therefore, there is a direct correlation between air temperature and the height of the boundary layer. Finally, the daily-averaged value of the boundary layer height for the whole 105 days of data was 2,067 km AGL. Also, the hourly-averaged height of the atmospheric boundary layer has been extracted and plotted for all four seasons. For these data, a good relationship is observed between air temperature and the height of the boundary layer, especially in winter, which is expected to improve accuracy if further measurements are made.

Geographical location of Iran in arid and semi-arid regions has strongly affected its food security, water resources and weather- and climate-related extreme events due to climate change. Global warming, on average, increases precipitation on Earth by increasing evaporation from the oceanic surface that enters the atmosphere, but the response of the West Asian region to global warming is generally a decrease in precipitation. Some studies have confirmed a decrease in the average precipitatiomn of Iran with an increase of precipitation in the south and southeast of the country. They also confirmed that the largest decrease in the precipitation of Iran occurs in the Zagros region.To this end, in this study, an overview of possible changes in precipitation trends in 43 stations of Iran is presented. To achieve the goals of the paper, the data of four Earth System Models from CMIP6 models, including MIROC6, FGOALS_g3, BCC-CSM2-MR and ACCESS-ESM1-5 were used. The precipitation output of the selected models has been statistically downscaled using CMHyd software for 43 synoptic stations over Iran. Observational period of 1985-2014 and the next three 25-year periods as the near future 2020-2026, the mid future 2075-2051 and the far future 2100-2076 were considered as study periods. Future changes in precipitation under three Shared Socio-economic Pathways (SSP) scenarios of SSP1-2.6, SSP2-4.5 and SSP5-8.5 were estimated. During the downscaling process, the Classius-Clapiron (CC) rate was applied to increase the trend of heavy rainfall. The air holding capacity is controlled by the Clausius-Clapiron relation, which is the water vapour holding capacity of the air at 7% /oC, the so-called Clausius–Clapeyron (CC) rate (bellow equation).Since the maximum temperature increase in the country was considered 10 degrees Celsius based on the worst possible scenario (SSP5-8.5), so all severe precipitation events were projected that were above the threshold of 130% compared to the normal period that were reduced to 130% of its normal amount.The results showed that future rainfall changes were not significant in about 78% of the stations. The increase and decrease of rainfall were significant in 19% and 3% of the stations, respectively. The largest increase in precipitation will occur in the south-southeast and the largest decrease will occur in the central Zagros area. The average rainfall of the country will increase by 0.4% annually (with a uncertainty range of 14%). On a seasonal scale, precipitation changes in spring, summer, autumn and winter were estimated to be + 15.2, -11, -6 and +3.5, respectively. Although a 15.2% increase in precipitation is projected for spring, the range of uncertainty with 81.9% indicates a lack of confidence in future precipitation estimates of spring season. The maximum uncertainty is related to the spring rains, which shows that the rains of this season are becoming more and more distrustful. Under warming conditions, more spring rainfall increases will occur. In summer, in the near future, summer rains tend to increase and then decrease with more level of Global Warming. After spring, most uncertainty is related to summer precipitation. The autumn precipitation tends to be less than normal with lower amount of uncertainty, although its decrease is not significant. The small range of the autumn rainfall uncertainty chart indicates the agreement of different model-scenarios in projection of autumn rainfall. Winter rainfall projection does not have significant uncertainty and almost all model-scenarios agree on the relatively low trend of rainfall increase. In this season, under higher amount of Global Warming, the precipitation increases and the amplitude of uncertainty also increases. Among the next three periods, the lowest and highest range of uncertainty in the annual changes of precipitation in the first and last decade of the present century will occur with the range of changes of 14 and 29.1%, respectively. The lowest and highest amplitude of uncertainty on the seasonal scale with 3.6 and 81.9 are estimated in winter (near future) and spring (far future), respectively. This situation indicates that in the future most precipitation fluctuations will occur in the spring.

Iran and its neighboring countries are located in an arid and semi-arid region, so they have been always affected by various climate disasters. Drought as one of the most important climatic hazards in the region, has attracted lots of researches. Drought could have negative impacts on different areas such as agriculture, water resources, industry, transportation and urbanization. Drought is basically a qualitative phenomenon that could have various definitions in different sectors. In a basic categorization, drought could be climatic, hydrological, agricultural and social droughts. Researches have always keen on defining appropriate indices for each kind of drought. Standardized Precipitation Index (SPI) and Standardized Precipitation and potential Evapotranspiration Index (SPEI) are good examples of climatic drought indices which have been widely used and popular between researchers.While monitoring the changes of these indices are a good way of drought analysis, effective drought management also needs accurate forecast of drought. A widely applied approach to forecast drought is through using the outputs of Coupled atmosphere- land- ocean General Circulation Models (CGCMs). CGCMs simulate the atmosphere-land-ocean system in a simplified way. Although, these models have been developed a lot in recent years, the forecasts are not well enough to be used in small regions. Thus, downscaling and post-processing methods have been widely applied to increase the resolution and accuracy of CGCM’s outputs, respectively.Downscaling methods could be in a way grouped into two categories: 1. Statistical, 2. Dynamical. Dynamical methods are based on physical theories and using parametrization and initial hypothesis. Statistical methods investigate the relationship between modeled variables and their corresponding observation in a period of time. This relationship is then applied to forecast that variable. Both mentioned methods have their own advantages and disadvantages. In recent years, combined dynamical and statistical approaches have been widely used. In this study, SPI and SPEI were forecasted using the outputs of CFS.v2. Firstly, the outputs of CFS.v2 in the reforecast period (1982-2010) were downscaled to 30km horizontal resolution using RegCM4. Due to limitation in computational power, the dynamical downscaling was limited to 1-month lead time. Then needed outputs for the rest of analysis were post-processed using Decision Tree (DT) and Support Vector Machine (SVM). The outputs of CFS.v2-RegCM4 system for precipitation and ERA5 were used as inputs and ideal outputs, respectively, of DT and SVM for calibration and validation. These post-processed outputs were then applied to calculate SPI. In order to calculate modeled SPEI, firstly, ET was computed using Hargreaves-Samani approach and the needed climatic modeled variables by CFS.v2-RegCM4 system. The reanalysis ERA5 data was utilized to calculate reanalysis ET. Then, in order to calibrate and validate DT and SVM models, modeled ET and reanalysis ET were applied as input and ideal output, respectively. Finally, the postprocessed precipitation and ET were employed to compute modeled SPEI. These modeled SPI and SPEI were spatially compared with SPI and SPEI computed using ERA5 data. The results showed that both DT and SVM improved the accuracy of precipitation and ET forecasted by CFS.v2-RegCM4 system. However, DT model performed much better than SVM. The ability of DT in post-processing precipitation was higher than it for ET. Comparing spatial pattern of modeled SPI and SPEI with reanalysis SPI and SPEI showed that DT model could replicate these values at an acceptable level of accuracy.

This study aimed to evaluate the performance of the WRF-Chem model in estimating the amount of NO2 and O3 in the Tehran region during recent summers (2019 to 2021). First, by investigating the NO2 and O3 concentrations (including hourly and daily averages) from the Air Quality Control stations in Tehran city, the days of maximum ozone concentration in summer were selected for which run the Weather Research and Forecasting Chemistry (WRF-Chem) model. The simulations were performed for 36-hours, and the first 12 hours were considered spin-up times. The simulations were conducted using two nests with 30 and 10 km resolutions, respectively, so that the second domain covered the Tehran region. Also, in the model settings, 35 levels were considered in the vertical direction, and the pressure at the highest level was 50 hPa. In these simulations, the Global Forecast System (GFS) data with a spatial resolution of 0.5 degrees and 6-hours time step was used as the initial and boundary conditions. Physical parameterization schemes that performed well in previous studies in simulating atmospheric pollutants dispersion were used in the model configuration. The Rapid Radiation Transfer Model (RRTM) and the Goddard Shortwave schemes were also used to simulate the long-and short-wave radiation, respectively. The Monin-Obukhov scheme was used to simulate surface layer fluxes and the Yonsei University (YSU) PBL scheme was also used to simulate boundary layer fluxes. The land surface fluxes were obtained from the NOAH model. In addition, the Grell and Devenyi ensemble scheme was used to parameterize moist convection and the WRF-Single-Moment-Microphysics 5-class scheme was used to parameterize microphysical processes. The RADM2 chemical mechanism was also used in this configuration. The PREP–CHEM–SRC emissions preprocessor (version 1.5) was used to produce anthropogenic, biogenic and biomass burning gridded emission over the user-specified simulation domain. The global emission data comes from RETRO and GOCART background emission data. Anthropogenic emissions of greenhouse gases and air pollutants including CO, NH3, NOx, SO2, NMVOC and CH4, were derived from the EDGAR_HTAP v5.0 Emissions Database with 0.1 horizontal resolution. Primary anthropogenic aerosol emissions of BC, OC and DMS from GOCART model databases were also used. Biogenic emissions were calculated using the MEGAN model. Also, 3BEM fire emissions which are prepared using PREP-CHEM-SRC are used. Evaluation of the results of WRF-Chem model simulations was performed by two methods of horizontal distribution and station evaluation. The evaluation results of the simulated horizontal distribution of NO2 and the one is taken from the OMI satellite data showed that considering the slight spatial displacement in the model results, the WRF-Chem model had good performance in simulating the maximum surface NO2 in all cases. Considering the spatial distribution of O3 in days with maximum ozone pollution, the model has simulated areas of maximum ozone in the Tehran. However, in most cases, no comment can be made on the accuracy of the simulated maximum areas because of the not precisely coinciding of the results maps with satellite transits. The station evaluation showed an overestimate of ozone concentration and a high underestimate of NO2 concentration by the WRF-Chem model.

Aerosols are solid or liquid particles suspended in the earth's atmosphere, which enter the earth's atmosphere from both natural and human sources. The wind blowing in the deserts, the evaporation of oceans and seas, the eruption of volcanoes, and the burning of forests and pastures are natural sources, and the burning of fossil fuels and the change of the earth's surface cover are human sources of their production. Aerosols can be classified into four types: dusty, marine, urban-industrial, and biomass burning particles. Due to the significant temporal and spatial changes of aerosols and for properly understand their climatic effects, we need to use long-term measurements of satellites and ground-based instruments. Measurements made from space and ground have allowed us to have a detailed view of the properties and effects of different types of atmospheric particles. Ground-based remote sensing is one of the powerful methods for determining the optical and physical properties of atmospheric aerosols.The sun-photometer (SPM) is a spectrometer that records the intensity of the sun's radiation usually in four wavelength channels of 440, 675, 870, and 1020 nm in two modes of measuring the sun and the sky with a limited viewing angle of 1.2 degrees during the day. Spectral aerosol optical depth, columnar water vapor, Angstrom exponent, single scattering albedo, polarized phase function, the real and imaginary refractive index of aerosols, and degree of linear polarization of sunlight are characteristics of atmospheric particles that are extracted from SPM measurements.There are different methods for classifying aerosols using data extracted from SPM measurements. One of the most common methods is to use aerosol optical depth data (a measure of the amount of atmospheric suspended particles) in terms of the Angstrom exponent (a qualitative measure of the dimensions of atmospheric particles). By combining other parameters obtained from the mode of the sun and the sky of the SPM, such as the aerosol optical depth, Angstrom exponent, particle size distribution, and refractive index, atmospheric aerosols can be classified. Our aim in this article is to investigate the ability of the degree of linear polarization parameter to classify the atmospheric particles. The degree of linear polarization measures the linear polarization of sunlight scattered by atmospheric particles (molecules and aerosols). For this purpose, the data of four sites of Banizoumbou, Beijing, El-Arenosillo, and Minsk, which have dusty, urban-industrial, marine, and biomass-dominant particles, respectively, were selected from the AERONET (AErosol RObotic NETwork) data.This paper uses three parameters of aerosol optical depth, Angstrom exponent, and degree of linear polarization extracted from SPM data. The results show that the maximum value of the degree of linear polarization (standard deviation) at the wavelength of 870 nm for dusty (Banizoumbou), urban-industrial (Beijing), marine (El-Arenosillo) and biomass (Minsk) aerosols are equal to 0.14 (0.05), 0.35 (0.10), 0.47 (0.08) and 0.37 (0.08) respectively. Therefore, the parameter of the degree of linear polarization is able to separate dusty, urban-industrial, and marine atmospheric particles from each other. However, biomass particles overlap a lot with urban-industrial aerosols and cannot be separated from each other.

In this research, three steps were taken to estimate the solar energy balance on the earth's surface. First, the amount of incident radiation on a tilted surface at the top of the atmosphere was calculated. Then, by using MODIS data, the transmittance coefficients of the atmosphere were estimated and the amount of direct radiation, diffuse radiation and global radiation in cloudless sky conditions were estimated. In the next step, based on the cloud transmittance coefficient, the amount of all sky radiation was estimated. Finally, by estimating the actual albedo of the earth's surface, the balance of solar radiation on the earth's surface was evaluated.The average top of atmosphere radiation in Iran is about 365 Watts per square meter. On a tilted surface, Iran receives 356 Watts per square meter of solar radiation. The difference in the angle of radiation on a tilted surface compared to the flat ground due to the slope of the ground and the difference in the duration of the radiation on a tilted surface compared to the flat ground due to the aspect of slope resulted a 2.5 percent reduction in the amount of radiation in Iran.In Iran, on a clear and sunny day about one percent of solar radiation is lost by air molecules not reaching the ground. The phenomenon of Rayleigh scattering also prevents about 9% of radiation from reaching the earth's surface. Therefore, about 10% of solar radiation is reduced due to atmospheric gases. The presence of aerosols, water vapor and ozone also affect the transparency of the atmosphere to solar radiation. The effect of these gases can be expressed by the transmission coefficient namely the aerosols transmittance coefficient which is low in desert areas of the country and on the coasts of Oman Sea and Persian Gulf and for Khuzestan Plain. In these areas, between 20 and 40 percent of the solar radiation is prevented from reaching the earth's surface by the aerosols. On the other hand, in the heights of Zagros and Alborz mountains and in the heights of Khorasan and in the north-west of Iran, aerosols do not play a significant role in reducing solar radiation. In Iran, the average reduction of solar radiation due to the presence of aerosols is about 17%.As expected, water vapor transmission is minimal at high altitudes, and about 10% of solar radiation is prevented from reaching the earth's surface due to atmospheric water vapor. On the shores of the Oman Sea, Caspian Sea, and Persian Gulf, the amount of attenuation due to atmospheric water vapor is about 14%. In Iran, the average reduction of solar radiation due to the presence of water vapor in the atmosphere is about 11%.The average transmittance of direct surface solar radiation in Iran is about 60%. In other words, the atmosphere prevents about 40% of direct sunlight from reaching the earth's surface. In mountainous areas the transmittance coefficient is the maximum and exceeds 70%. In the southern banks and eastern and central regions of Iran, due to the presence of aerosols and water vapor, the figure is less than 60%. The amount of mean direct radiation in Iran is about 213 Watts per square meter. Diffuse radiation is a small part of the total radiation. The average transmittance of diffuse radiation in Iran is about 10%. Aerosols play an important role in scattering solar radiation. The amount of mean diffuse radiation that reaches the earth's surface in Iran is about 35 Watts per square meter.This study shows that the global radiation in Iran is 248 Watts per square meter. The average transmittance coefficient of global radiation is 70% and follows the configuration of topography and distance from the sea. Average cloudiness of Iran is about 26% and the average ratio of actual to possible sunshine hours is about 72%. On the shores of the Caspian Sea, the cloudiness exceeds 60%. The average cloud transmittance coefficient in Iran is about 83%. In Iran, clouds contribute about 17% in the reduction of radiation. On a cloudy day, the mean amount of solar radiation that passes through the atmosphere and reaches the surface of the earth on a tilted surface is 205 Watts per square meter. The average albedo of Iran is about 21%. Nearly 80% of the solar radiation that reaches the earth's surface is absorbed by the surface. The amount of net annual solar radiation on the earth's surface in Iran varies between 80 and 220 Watts per square meter.

One of the effects of drought is the loss of vegetation, drying of wetlands, and desertification, which in turn can create an appropriate substrate for the dust and increase the concentration of aerosol. In this research, with the aim of understanding the spatial and temporal changes of dust storms in Golestan province over a period of 20 years (2000 to 2020), AOD (Aerosol Optical Depth) data from MODIS were used. So, dust detection when the AOD value was greater than 0.5, was considered a dust event. Also, the monthly horizontal visibility data of 10 synoptic stations with appropriate distribution in the province were used to validate the AOD product. Monthly horizontal visibility data is obtained by averaging three-hour horizontal visibility at 10 synoptic stations in Golestan province. With the aim of investigating the annual and seasonal distribution of dust in the Golestan province, a diagram of the number of days with dust in the study time was obtained. Then, the map of the spatial changes of dust was also presented. To investigate the intensity of influence of Golestan province from the deserts of Turkmenistan, in addition to tracking the dust transfer path using the output of the HYSPLIT model for a selected event from each season, synoptic analysis was also performed using ERA5 data. Also, the trend of dust changes in the Golestan province and deserts of Turkmenistan was investigated. The validation results with a correlation of 0.66 show the ability of AOD data to present the dust distribution in Golestan province. According to the obtained results, in fact, there are two time periods with low occurrences of a dust storm (2000 to 2007 (except 2003) and 2016 to 2020) and one-time period with high occurrences of a dust storm (2008 to 2015). The highest number of days with dust is in the summer season and the lowest is in the winter season. In terms of spatial distribution, the northern regions of the province, i.e., Gonbad-e kavus, Aqqala, and Gomishan, are more involved in dust phenomena than other parts of the province. During the statistical period, the highest AOD value occurred on 7th September, 2020 with a value of 4.1. Tracing the dust transfer path using the output of the HYSPLIT model in the region showed that the main source of dust in the province was from Turkmenistan and especially the deserts of the Balkan province. Tracking the dust transfer route using the output of the HYSPLIT model in the region showed that the main source of dust in the province is mainly from the desert areas of Turkmenistan. The synoptic analysis of the selected day from each season generally indicates a strong pressure gradient in the region of Turkmenistan and north of Iran in such a way that the northerly winds in Golestan and the east of the Caspian Sea are significant, which causes the soil to rise, and it is transported to Golestan province. Also, on the selected days of spring, summer, and autumn, there is a ridge pattern and stable atmospheric conditions in the middle level of the troposphere, leading to dust persistence. In the winter, due to the presence trough of low in the middle of the troposphere and the occurrence of precipitation, the intensity of the dust is weaker. Also Comparing the annual and seasonal distribution of dust storms in the Golestan province and Turkmenistan shows the simultaneity of events.